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Indirect  Dominant  Mode  Rejection:  A  Solution  to 
Low  Sample  Support  Beamforming 

Ernesto  L.  Santos,  Michael  D.  Zoltowski,  Fellow,  IEEE ,  and  Muralidhar  Rangaswamy,  Fellow,  IEEE 


Abstract — Under  conditions  of  low  sample  support,  a  low-rank 
solution  of  the  minimum  variance  distortionless  response 
(MVDR)  equations  can  yield  a  higher  output  signal-to-in- 
terference-plus-noise  ratio  (SINK)  than  the  full-rank  MVDR 
beamformer.  In  this  paper,  we  investigate  several  low-rank  beam¬ 
forming  techniques,  and  we  also  propose  a  new  beamformer  that 
we  refer  to  as  the  indirect  dominant  mode  rejection  (IDMR).  We 
analyze  the  degradation  in  the  output  SINR  caused  by  residual 
cross  correlations  embedded  in  the  sampled  covariance  matrix 
due  to  low  sample  support.  The  IDMR  beamformer  is  based  on  a 
parametric  estimate  of  the  covariance  matrix,  in  which  any  cross 
correlation  is  canceled  out.  Simulations  reveal  that  the  IDMR 
beamformer  yields  a  dramatic  improvement  in  output  SINR 
relative  to  the  conjugate  gradient  (CG),  principal  component 
inverse  (PCI),  and  dominant  mode  rejection  (DMR)  beamformers. 
In  our  investigation  of  the  low-rank  CG  beamformer,  we  address 
the  issue  of  whether  the  unity  gain  constraint  in  the  look  direction 
should  be  enforced  a  priori  via  the  use  of  a  blocking  matrix  or 
effected  a  posteriori  through  simple  scaling  of  the  beamforming 
vector.  Remarkably,  it  is  proven  that  the  two  methods  yield  exactly 
the  same  low-rank  beamformer  at  each  and  every  rank. 

Index  Terms — Blocking  matrix,  conjugate  gradient  beamformer, 
dominant  mode  rejection  (DMR)  beamformer,  indirect  dominant 
mode  rejection  (IDMR)  beamformer,  low  sample  support  beam¬ 
forming,  minimum  variance  distortionless  (MVDR)  beamformer, 
robust  adaptive  beamforming. 


I.  Introduction 

ARRAYS  of  sensors  have  been  widely  used  in  wireless 
communications  in  a  variety  of  applications  such  as 
sonar,  radar,  astronomy,  medical  imaging,  sound  processing, 
and  seismic  exploration.  In  this  paper,  we  investigate  several 
low-rank  beamforming  techniques,  and  we  also  propose  a  new 
beamforming  technique  that  we  refer  to  as  indirect  dominant 
mode  rejection  (IDMR)  [1].  To  maximize  performance  in  a 
nonstationary  scenario,  adaptive  beamforming  (ABF)  is  nec¬ 
essary.  In  rapidly  changing  scenarios,  only  a  limited  number 
of  snapshots  is  available  to  estimate  the  covariance  matrix; 
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under  these  circumstances,  the  low-rank  minimum  variance 
distortionless  response  (MVDR)  equations  yield  faster  con¬ 
vergence  in  terms  of  the  signal-to-interference-plus-noise  ratio 
(SINR)  performance  than  the  full-rank  MVDR  beamformer 
[2]— [5] .  Owsley  [6]  introduced  dominant  mode  rejection 
(DMR)  low-rank  ABF.  Kirsteins  and  Tufts  [7]  introduced 
principal  components  inverse  (PCI)  low-rank  ABF.  Goldstein 
et  al.  [8]  introduced  the  multistage  Wiener  filter  (MWF),  which 
was  proven  to  be  equivalent  to  the  conjugate  gradients  (CG) 
ABF  [9].  The  performance  of  the  aforementioned  beamformers 
depends  on  the  rank  in  which  they  operate.  However,  there  is 
no  robust  algorithm  to  select  the  optimal  rank  because  the  true 
interference-plus-noise  covariance  matrix  is  not  known.  It  is 
well  known  that  the  optimum  rank  is  a  function  of  the  number  of 
signal  sources  (dominant  modes)  and  of  the  number  of  samples 
used  to  estimate  the  covariance  matrix.  Through  simulations  we 
also  show  that  the  optimum  rank  varies  with  the  look  direction 
[10],  i.e.,  the  optimal  rank  of  operation  is  not  necessarily  the 
same  for  all  the  different  signals.  In  the  literature  [6],  [11], 
the  DMR/PCI  method  is  applied  using  all  of  the  dominant 
eigenvectors  of  the  covariance  matrix.  The  number  of  dominant 
eigenvectors  does  not  vary  with  look  direction,  and  therefore 
does  not  affect  all  look  directions  with  an  optimal  SINR.  In 
our  simulations,  we  show  that  the  output  SINR  can  undergo 
substantial  variations,  depending  on  the  rank  that  is  chosen. 
Consequently,  the  output  SINR  can  be  seriously  degraded  if 
the  beamformer  operates  at  an  inadequate  rank.  To  solve  signal 
mismatch  and  low  sample  support  problems,  Vorobyov  et  al. 
proposed  in  [12]  a  robust  adaptive  beamforming  algorithm 
using  worst  case  performance  optimization.  Also  dealing  with 
signal  mismatch  and  low  sample  support  problems,  Mestre  and 
Lagunas  in  [13]  propose  an  optimal  diagonal  loading  technique. 
In  this  paper,  we  contribute  to  the  ABF  field  by  proposing  the 
IDMR  beamformer,  which  operates  at  its  maximum  rank  with 
an  output  SINR  very  close  to  the  maximum  over  all  ranks. 
Also,  the  output  SINR  is  significantly  higher  than  the  SINR 
obtained  with  existing  ABF  techniques. 

The  MVDR  beamformer  maximizes  the  output  SINR  only 
when  the  desired  signal  is  uncorrelated  to  the  interference. 
Due  to  finite  sample  averaging,  residual  correlations  between 
sources  are  present  in  the  sample  covariance  matrix,  causing  a 
degradation  in  the  performance  of  MVDR-based  beamformers. 
Even  when  all  incident  signals  are  truly  uncorrelated,  the 
sample  signal- source  covariance  matrix  intrinsically  embedded 
in  the  sample  spatial  covariance  matrix  has  substantial  off-di¬ 
agonal  terms.  This  is  especially  true  under  conditions  of  low 
sample  support,  thereby  giving  the  “appearance”  of  correlation 
between  signals.  This  problem  is  circumvented  in  the  IDMR 
algorithm  proposed  herein  by  forcing  the  correlation  between 
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signals  to  be  zero  in  the  parametrically  constructed  spatial  co- 
variance  matrix  estimate,  which  is  then  used  in  the  computation 
of  the  MVDR  beamformer.  Through  simulations,  we  show 
that  under  low  sample  support,  the  IDMR  beamformer  yields 
a  substantially  higher  SINR  than  either  the  full-rank  MVDR 
beamformer  or  low-rank  MVDR  beamformer s  such  as  CG  and 
PCI/DMR,  even  when  these  low-rank  beamformers  operate 
at  the  optimal  rank.  We  also  show  through  simulations  that 
when  signals  are  truly  correlated,  the  IDMR  beamformer  yields 
substantially  higher  SINR  than  the  low-rank  CG  and  PCI/DMR 
beamformers.  Even  though  the  IDMR  beamformer  was  orig¬ 
inally  designed  as  a  tool  to  circumvent  residual  signal- source 
correlation  due  to  low  sample  support,  it  is  also  effective  in  the 
case  where  the  sources  are  truly  correlated. 

The  relationship  between  steering-independent  ABF 
and  steering-dependent  ABF  is  investigated  for  the  case  of 
CG-based  low-rank  beamforming.  The  scenario  assumes  the 
formation  of  multiple  adaptive  beams,  each  pointed  to  a  dif¬ 
ferent  look  direction.  In  steering-dependent  ABF,  a  generalized 
sidelobe  canceler  (GSC)  [14]  is  formed  for  each  look  direction. 
Mathematically,  the  GSC  serves  to  convert  the  constrained 
MVDR  optimization  problem  to  an  unconstrained  optimization 
problem,  thereby  enforcing  a  priori  the  unity  gain  constraint 
in  the  look  direction.  In  contrast,  in  steering-independent 
ABF,  an  unsealed  version  of  the  Wiener-Hopf  equations  is 
solved,  and  the  unity  gain  constraint  is  enforced  a  posteriori 
through  simple  scaling  of  the  resulting  ABF  weight  vector. 
Implementation  of  a  GSC  for  each  look  direction  requires  the 
construction  and  application  to  the  data  of  a  blocking  matrix  for 
each  look  direction.  The  attendant  computational  complexity  is 
quite  substantial.  In  this  paper,  we  prove  a  very  important  and 
somewhat  surprising  result:  the  low-rank  beamformer  obtained 
with  steering-dependent  (SD)  conjugate  gradients  (SD-CG)  is 
exactly  the  same  as  the  low-rank  beamformer  obtained  with 
steering-independent  (SI)  conjugate  gradients  (SI-CG)  at  any 
rank.  Thus,  the  performance  of  SD-CG  can  be  obtained  without 
having  to  form  blocking  matrices  for  each  look  direction. 

In  Section  II,  we  present  the  signal  model.  In  Section  III, 
we  present  some  background  on  MVDR  beamforming  and 
provide  an  overview  of  the  SI-CG  and  SD-CG  beamformers. 
In  Section  IV,  we  prove  the  equivalence  of  the  SI-CG  and 
SD-CG  beamformers.  In  Section  V,  we  show  the  deterioration 
to  MVDR-based  beamformers  caused  by  residual  correlations, 
and  we  introduce  the  IDMR  beamformer.  In  Section  VI, 
we  present  the  simulations  where  we  compare  the  IDMR 
beamformer  to  the  existing  beamformers.  In  Section  VII,  we 
summarize  the  conclusions. 

II.  Signal  Model 


input.  We  will  also  assume  additive  noise  at  each  sensor.  Fet 
x(t)  =  [%i(t),X2(t),  •  •  •  where  xt (t)  is  the  signal 

received  at  the  ith  array  element.  (.)T  denotes  transpose.  Given 
the  above-assumed  scenario,  x(t)  can  be  expressed  in  the 
following  fashion: 

x(t)  =  As  (t)  +  n(t).  (1) 

A  =  [a(#i),  a(#2),  •  •  • ,  a(0d)]  is  called  the  signal  direction  ma¬ 
trix  (SDM),  where  Aik  =  ai(0k)  is  the  response  of  the  it h  array 
element  relative  to  that  of  the  first  element  when  a  single  signal 
arrives  at  Ok ,  s (t)  =  ,  Sd(t)]T,  where  Sk(t)  is  the 

signal  associated  with  the  &th  source  as  received  at  the  first  array 
element,  and  n (t)  =  [ni(Z),  712(f),  •  •  • ,  nrn{t)]T ,  where  rii(t)  is 
the  additive  noise  present  at  the  it h  array  element. 

The  output  of  the  beamformer  is  given  by  y(t)  =  w Hx(t), 
where  w  =  [uq,  W2, . . .  ,  wm]T  is  the  beamformer  weight 
vector.  (.)H  denotes  Hermitian  (conjugate)  transpose. 


III.  MVDR  Beamforming 

After  beamforming,  the  output  SINR  can  be  expressed  as 
[12],  [15] 


SINR  = 


£j|wHaiSi(X)|2j 

£{|wh[x(T)  -aisi(i)]|2} 


W/fR;  1,,  W 


where 


Ri+n  =  £{[x(0  -aisi(^)][x(t)  -aisi(t)]ff  j  (3) 

is  the  interference-plus-noise  covariance  matrix,  and 
<j\  —  E{\si(t)\2}  is  the  power  of  the  desired  signal.  Maxi¬ 
mizing  (2)  is  equivalent  to  minimizing  1/SINR,  as  follows: 

wffRi+nw 

111111  21 - F7 - io-  (4) 

w  o-f\wna.i\z 

Observing  that  the  cost  function  above  is  unaltered  when  w 
is  scaled,  we  can  impose,  without  a  loss  of  generality,  the  con¬ 
straint  w^ai  =  1,  which  can  always  be  achieved  by  scaling  w. 
Hence,  the  maximization  of  the  output  SINR  is  equivalent  to 

min  w^R7;+nw 

w 

subject  to  :  w^ai  =  1.  (5) 

The  solution  to  (5)  is  the  well-known  MVDR  beamformer, 
given  by  solving 


Consider  an  array  of  m  sensors  receiving  signals  from  d 
sources  of  emission  at  directions  Ok,  k  =  1, . . . ,  d  with  re¬ 
spect  to  the  broadside  of  the  array.  The  angles  of  arrival  are 
enumerated  such  that  <9i  is  the  angle  of  the  desired  source. 
We  will  assume  that  the  desired  source  is  narrowband  and 
that  the  narrowband  filtering  about  the  center  frequency  of  the 
desired  source  fQ  occurs  at  the  front  end  of  the  receiver  such 
that  the  (d  —  1)  interfering  signals  are  narrowband  and  colo¬ 
cated  in  frequency  with  the  desired  signal  at  the  beamformer 


Ri+nW  =  Aai  (6) 

where  A  serves  to  satisfy  the  unity  gain  constraint  in  (5). 

In  many  practical  applications,  the  ideal  signal-free  covari¬ 
ance  matrix  R;+n  is  unavailable,  and  a  sample  estimate  from 
N  snapshots 

1  n 

R  =  -  F]x(7)xh(7)  (7) 

n  t= 1 
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is  often  used  instead  of  R^+n.  It  is  observed  that  if  the  desired 
signal  is  uncorrelated  to  the  interference  and  to  the  array  noise, 
then 


min  w^Rw 

w 

subject  to  :  ai  =  1  (8) 

where 

R  =  E{x(t)*H(t)}  (9) 

is  equivalent  to  the  minimization  problem  of  (5).  Observe  that 
under  the  constraint  w^ai  =  1 


min  w^Ri+nw  =  minw H  E 

w  w 

x  |  [As(i)  -  aisi(i)  +  n(i)] 

x  [As(i)  -  aisi(i)  +  n(f)]H  j  w 
=  min  {wffRw  —  a\  |wffai|2} 

W  1  J 

=  min  wffRw.  (10) 

w 

Hence,  when  the  desired  signal  is  uncorrelated  to  the  interfer¬ 
ence  and  to  the  array  noise,  minimizing  w^R^+nw  subject 
to  w^ai  =  1  is  the  same  as  minimizing  w^Rw  subject  to 
w^ai  =  1;  and,  therefore,  the  MVDR  beamformer  with  ei¬ 
ther  R^_|_n  or  R  maximizes  the  output  SINR.  It  is  important  to 
note  that  this  only  happens  when  the  true  covariance  matrix  is 
known  and  the  desired  signal  is  uncorrelated  to  the  interference 
and  to  the  array  noise.  Even  when  sources  are  truly  uncorre¬ 
lated,  residual  correlations  arise  due  to  finite  sample  averaging, 
causing  the  output  SINR  to  be  significantly  degraded,  since  min¬ 
imizing  wfRw  is  not  tantamount  to  minimizing  w^R^+nw. 
This  motivates  us  to  introduce  IDMR  beamforming,  which  is 
based  on  a  parametric  estimate  of  the  correlation  matrix  where 
the  residual  correlations  are  canceled  out.  The  IDMR  beam- 
former  is  presented  in  Section  V. 

Low-rank  MVDR  steering-independent  beamformers  arise 
when  a  a  low-rank  solution  for  (6)  is  used.  SI-CG  consists 
of  using  the  CG  algorithm  [16]  to  solve  (6);  this  solution  is 
obtained  by  minimizing  the  quadratic  function 

/( w)  =  ^w/fRw  —  w^ai .  (11) 

Steering-dependent  ABF  consists  of  enforcing  the  unity  gain 
constraint  a  priori  through  the  use  of  a  blocking  matrix  B  and 
thus  transforming  the  constrained  optimization  problem  of  (5) 
in  an  unconstrained  optimization  problem.  The  columns  of  B 
form  an  orthonormal  basis  for  the  orthogonal  complement  of 
the  steering  vector  ai .  While  B  is  not  unique,  it  must  satisfy 

K(  B)  =3^ 

B^B=I  (12) 

where  the  operator  1Z  denotes  range  of  a  matrix,  and  I  denotes 
the  identity  matrix.  There  are  several  forms  to  find  a  blocking 


matrix  B;  computation  through  a  Householder  transformation 
is  shown  in  the  Appendix.  Expressing  w  as 

w  =  —  (ai  +  Bu)  (13) 

m 

guarantees  that  the  unity  gain  constraint  is  satisfied  a  priori. 

Substituting  (13)  into  (5)  leads  to  the  unconstrained  optimiza¬ 
tion  problem 

min  -Ar(ai  +  Bu)^R(aj  +  Bu).  (14) 

u  mz 

Taking  the  gradient  of  (14)  with  respect  to  u  dictates  that  the 
optimal  value  of  u  is  the  solution  to 

(B^RB)u  =  -B^Rai.  (15) 

Low-rank-MVDR- steering-dependent  adaptive  beamforming 
methods  consist  of  using  a  low-rank  solution  for  (15)  and 
substituting  it  in  (13).  When  the  full-rank  solution  is  taken,  both 
solutions  of  (6)  and  (13),  SI  and  SD  beamforming,  respectively, 
yield  the  same  weight  vector  w.  When  a  low-rank  solution  is 
used,  the  low-rank  SI  and  SD  weight  vectors  are  not  guaranteed 
to  be  equal. 

IV.  Equivalence  Between  SI-CG  and  SD-CG 

In  the  SI-CG  beamformer,  the  CG  algorithm  is  used  to  solve 
(6),  note  that  (see  [17]  and  [18])  at  step  r,  the  output  vector  of  the 
CG  algorithm  w minimizes  the  R-norm  of  the  error,  which 
is  given  by 

esi  (w'M)  =  (w'W  -  R^ai^R  (w'W  -  R-V) 

(16) 

subjected  to  being  constrained  within  the  Krylov  subspace 

Kg,  defined  by  the  column  space  of  below: 

=[a1iRa1i  ...  iR’-V].  (17) 

It  follows  from  [9]  that  after  r  steps  of  CG,  the  solution  is  con¬ 
strained  to  lie  in  a  r-dimensional  Krylov  subspace.  Thus,  it  is 
referred  to  as  a  rank  r  solution. 

For  the  purposes  of  comparing  the  performance  of  SD-CG 
relative  to  SI-CG,  we  desire  that  the  rank-one  solution  for  both 
methods  be  the  scaled  steering  vector  (l/m)ai,  which  is  tanta¬ 
mount  to  the  CBF.  Thus,  for  SD-CG,  (r  —  1)  steps  of  CG  are 
taken  in  solving  (15)  to  obtain  the  solution  denoted  u(r_1) .  Sub¬ 
stituting  u(r_1)  for  u  into  (13)  yields  the  rank  r  SD-CG  beam- 
former  denoted  w^J . 

Similarly  to  w'(r) ,  u(r_1)  is  restricted  to  lie  in  the  Krylov  sub¬ 
space  JCu~  ^  defined  by  the  range  space  of  the  matrix  K.u~1\ 
where  is  defined  as 


Kir ) 


B^Rai  :  (B^RB)B^Rai  : 


:(BHRB)r  1BHRai 
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(T) 

In  turn,  (13)  and  the  equation  above  dictate  that  wKsdJ  is  restricted 
to  lie  in  the  range  space  of  the  matrix 


K 


(r)  _ 

sd 


(18) 


Simple  algebraic  manipulation  leads  to  the  following  expression 
for  the  Krylov  subspace  represented  by  the  range  space  of 


K 


( r )  _ 
sd 


ai  :  (BBHR)ai  : 


:(BB^R)r  'a. 


(19) 


The  following  theorem  states  that  the  SD-CG  and  SI-CG 
Krylov  subspaces  of  dimension  r  are  the  same. 

Theorem  1: 

spanjai ,  Rai , . . . ,  Rrai }  =  spanjai ,  (BBHR)ai , 

. . . ,  (BB^R)rai},  r  —  1,2,....  rri. 


Proof  is  by  induction.  The  basis  case  is  proven  in  step  1,  and 
the  induction  step  is  proven  in  step  2  below. 

1)  Basis  case: 

span{ai,Rai}  =  spanjai,  BB^Rai}. 


Proof. 

i)  [ai:B]  is  an  orthonormal  basis  for  IZ171  since  B  is 
the  orthogonal  complement  of  ai  by  construction, 
and  the  columns  of  B  are  orthonormal.  BBH  is  a 
projection  matrix. 


spanjai,  Cai, . . . ,  Crai,  Cr+1ai} 

=  spanjai,  Rai, . . . ,  Rrai,  Cr+1ai}, 
by  induction  hypothesis 
=  spanjai ,  BB^  Rai , . . . ,  BBHRrai , 

r 

^  G^BBjffRt+1ai},  by  i)  and  ii) 

7-0 

=  spanjai,  BB^Rai, 

. . . ,  BB^Rrai,  BB^Rr+1ai} 

=  spanjai,  Rai, . . . ,  Rrai,  Rr+1ai},  by  i). 

Proof:  ■ 

Thus,  the  rank  r  SI-CG  and  SD-CG  beamformers  are  con¬ 
strained  to  lie  in  the  same  subspaces. 


(r) 

A.  Closed-Form  Expressions  for  wysiJ 

Since  we  have  proved  that  the  rank  r  beamformers  for  SI-CG 
and  SD-CG  lie  in  the  same  subspace,  the  unsealed  (prior  to 
scaling  to  satisfy  the  unity  gain  constraint)  SI-CG  beamformer 
Wgf  may  be  expressed  as 


(20) 


.*.  3/?|Rai  =  /3ai  +  BB^Rai. 

ii)  span{ai,Rai}  =  spanjai ,/iai  +  BB^Rai}= 
spanjai,  BB^Rai}  ■ 

2)  Inductive  step: 

Define  C  =  BB^R. 

If  spanjai ,  Rai , . . . ,  Rrai }  =  spanjai ,  Cai , . . . , 
Crai},  then  spanjai,  Rai, . . . ,  Rr+1ai}  = 
spanjai,  Cai,  . . . ,  Cr+1ai} 

Proof. 

i) 

V  j  >  0,3/3j|RJai 

=  6j8l\  +  BBHRJai, 

Vfc  >  0,  spanjai,  Rai, . . . ,  Rfcai} 

=  spanjai, /?iai  +  BBHRai, 

. . . ,  /3fcai  +  BB^R^ai} 

=  spanjai,  BB^Rai, . . . ,  BB^R^ai}. 

ii) 

Cr+1ai 

=  C(Crai) 

r 

=  c(£  a/R'ai),  by  induction  hypothesis 
7-0 

7-0 


where  y(r)  and  are  r  x  1  vectors  whose  elements  com¬ 
bine  the  columns  of  Ks^  and  Ks^,  respectively.  At  this  point, 
we  seek  a  closed-form  expression  for  x(r) ,  with  given  by 
(19).  With  this,  the  rank  r  SI-CG  beamformer  is  expressed  in 
the  Krylov  basis  of  the  SD-CG  beamformer.  Later,  we  find  a 
closed-form  expression  for  the  rank  r  SD-CG  beamformer  in  the 
same  basis.  We  can  then  easily  compare  the  SI-CG  and  SD-CG 
beamformers. 

The  value  of  x^r)  is  obtained  by  substituting  = 

K^XW  into  the  objective  function  in  (16)  and  noting  that  x(r) 
minimizes  the  function.  Subsequently,  taking  the  gradient  and 
setting  it  to  zero  yields 


K 


(r)ffRK 


sd 


KS'%  =  0 


=  (KfRK^J^KfaL 


(21) 


Substituting  the  expression  for  in  (18),  the  inverse  of 
K^RK^  may  be  expressed  as 


(k^RK«)_1 

afRai  afRBK^"1) 

“  _Klr_1)ffBHRai  Klr_1)ffBHRBKlr-1) 

In  addition,  it  follows  that 


Kifa,  =  [1  i  0]» 


(23) 
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Substituting  (22)  and  (23)  into  (21)  dictates  that  only  the  first 
column  of  the  inverse  matrix  in  (22)  is  needed.  Invoking  the 
block  matrix  inversion  lemma,  may  be  computed  as 

X(  }  ~~  7  -  (Kir_1)HBffRBK(lr_1))_1Kir_1)ffBRai. 

(24) 

where  7  is  a  multiplicative  scalar  that  will  be  accounted  for 
when  the  unity  gain  constraint  is  satisfied.  Substituting  (24)  into 
(20),  with  given  by  (18),  yields 

wI7  =  7  (aa  -  BK^1)  (k^-^B^RBKJT1))-1 

xKf-^BRai)  .  (25) 

Since  the  columns  of  B  are  orthogonal  to  ai  by  construction, 
scaling  to  satisfy  the  unity  gain  constraint  simply  implies 

=  2  [ai  -  BKlT1)  (Kr^B-RBKr1))-1 

xKff_1)ffBRai  .  (26) 

B.  Closed-Form  Expressions  for 

According  to  (13),  the  rank  r  SD-CG  beamformer  may  be 
expressed  as 

w 7  =  1  [a,  -  Bu^-1)’ 

sa  m  L  J 

=  -  [ai  +  BK',;'-|iz*,'-|)l  .  (27) 

ml  J 

Since  u(r_1)  is  determined  by  solving  (15)  via  the  CG  algo¬ 
rithm,  it  is  the  vector  that  minimizes  the  objective  function 

esd  (u(r-1))  =  (u<r-1>  +  (BffRB)_1BffRai)ff 

xB^RB  (ul>  1:1  +  (BffRB)_1BffRai)  (28) 

with  u(r_1)  constrained  to  lie  in  the  range  space  of  the  ma¬ 
trix  K^-1)  in  (18).  Substituting  Ku~^z^r~^  foru(r_1)  above 
dictates  that  z(r_1)  is  the  vector  that  minimizes 

z¥^HK^HBHRBKf^h^r-^ 

+  af  RBK^^z^^ 


Equations  (26)  and  (30)  are  observed  to  be  identical,  thereby 
proving  that  the  rank  r  SI-CG  beamformer  is  equal  to  the  rank  r 
SD-CG  beamformer.  Brennan,  Mallett,  and  Reed  have  shown  in 
[19]  that  the  full-rank  MVDR  steering-dependent  beamformer 
implemented  through  a  sidelobe  cancellation  [14]  and  the  full- 
rank  MVDR  steering  independent  beamformer  converge  to  the 
same  beamformer.  We  showed  that  when  implemented  through 
CG  the  low-rank  versions  of  these  beamformer  also  yield  the 
same  beamformer  at  every  rank. 

V.  Indirect  DMR 

In  Section  II,  we  showed  that  when  the  desired  signal  is  un¬ 
correlated  to  the  interference,  the  MVDR  beamformer  maxi¬ 
mizes  the  output  SINR.  This  is  true  even  when  the  desired  signal 
is  present  in  the  covariance  matrix,  as  long  as  the  true  correla¬ 
tion  matrix  is  available.  If  the  desired  signal  is  correlated  with 
the  interference,  then  the  MVDR  beamformer  no  longer  max¬ 
imizes  the  output  SINR.  Even  when  the  desired  signal  and  the 
interference  are  not  correlated,  the  covariance  matrix  estimated 
from  a  limited  number  of  snapshots  contains  residual  correla¬ 
tion.  Under  these  circumstances,  the  MVDR  beamformer  does 
not  maximize  the  output  SINR  [20]. 

As  an  example,  we  consider  the  case  in  which  two  uncorre¬ 
lated  signals  are  arriving  at  the  array.  To  facilitate  the  analysis, 
we  do  not  include  the  effects  of  the  array  noise.  The  signal  ar¬ 
riving  at  the  array  is  given  by 

x(t)  =  Si(i)ai  +  S2(t)  a2  (31) 

where  $%(t)  and  s2(t)  are  independent  complex  Gaussian 
random  processes.  The  sample  covariance  matrix  from  n  sam¬ 
ples  is  given  by 

1  n 

=  [»i  a2] 

x*  (32) 

n-lY^iS2{f)st{t)  n~^U\s2{t)\2  J  I  a2  J 

The  nondiagonal  terms  in  the  signal- source  covariance  matrix 
(SSCM),  which  is  represented  by  the  matrix  in  the  middle,  con¬ 
verge  to  zero  as  n  — ►  oo.  This  is  proved  by  observing  that  the 
mean  of  these  nondiagonal  terms  is  zero 

£(jEs2(fK(4=0  (33) 


where  we  have  omitted  an  inconsequential  additive  constant. 
Taking  the  gradient  with  respect  to  n(r_1)  and  setting  it  to  zero 
yields 

z(r_1>  =  -  ^K[’’_1*iJBffRBK[r"1^"1  K^r_1^BffRai. 

(29) 

In  turn,  substituting  z(r_1)  back  into  (27)  yields 
W {;J  =  2  [ai  -  BK71)  (Kf^B^RBKr1))-1 

xK^-'^B^Ra,]  .  (30) 


and  the  variance  goes  to  zero  as  n  — >  oo 


e  {  t  it  it 


where  o\  —  E?{|si(£)|2}  and  cj\  —  E{\si(t)\2}. 

Since  the  variance  of  the  nondiagonal  terms  go  to  zero  rela¬ 
tively  slowly  (proportional  to  1  /n),  under  low  sample  support, 
these  terms  can  significantly  deteriorate  the  output  SINR  ob¬ 
tained  with  MVDR  beamforming.  This  motivates  us  to  propose 
the  IDMR  beamformer,  which  is  based  on  a  parametric  estimate 


3288 


IEEE  TRANSACTIONS  ON  SIGNAL  PROCESSING,  VOL.  55,  NO.  7,  JULY  2007 


TABLE  I 

Outline  of  the  IDMR  Algorithm 


i 

Estimate  the  directions  of  the  dominant  sources  and  form  A. 

a 

Estimate  . 

Hi 

Compute  —  diag  (At  (R  —  A^H^j . 

iv 

Parametrically  form  R idmr  —  AP^A^  +  I. 

V 

Compute  Widmr  =  ( R,>W  al )  “ 1  RUmr  al 

of  the  covariance  matrix  in  which  the  embedded  residual  corre¬ 
lations  are  canceled  out.  In  the  derivation  of  the  IDMR  beam- 
former,  we  assume  that  the  array  manifold  is  known.  In  some 
practical  applications,  there  is  an  error  in  the  estimate  of  the 
array  manifold.  IDMR  beamforming  in  the  presence  of  array 
manifold  mismatch  will  be  considered  in  future  work. 

Substituting  (1)  into  (9)  and  assuming  that  s (t)  and  n(t)  are 
uncorrelated  and  E{n(t)nH (t)}  =  cr^I,  where  is  the  power 
of  the  noise  at  each  array  element,  then  the  true  covariance  ma¬ 
trix  R  can  be  given  by 

R  =  APAh  +  a2n  I  (35) 

where  the  matrix  P  is  the  signal  source  covariance  matrix  and 
is  given  by 


P  =  E  {s(t)sH (t)}  .  (36) 


In  the  IDMR  beamformer,  we  propose  to  use 

Pd  =  diag  (At  (R  -  a2nl)  A*H)  (37) 


as  an  estimate  of  the  diagonalized  signal  source  covariance  ma¬ 
trix,  where  A”*"  =  (A^ A)  AH  is  the  pseudoinverse  of  A. 
A  is  formed  from  the  estimated  signal  directions,  given  the 
known  form  of  the  array  manifold,  and  is  an  estimate  of  the 
noise  variance.  The  nondiagonal  elements  of  P,  arising  from 
residual  correlations  among  sources  due  to  finite  sample  aver¬ 
aging  (or  even  due  to  true  correlation  among  sources)  are  dis¬ 
carded,  giving  rise  to  a  diagonalized  signal- source  covariance 
matrix  estimate  denoted  TV  Note  from  Section  II  that  the  nondi¬ 
agonal  elements  cause  the  MVDR  beamformer  to  not  maximize 
the  output  SINR.  Thus,  even  when  sources  are  truly  correlated, 
discarding  these  terms  allows  the  MVDR  beamformer  to  maxi¬ 
mize  the  output  SINR. 

To  compute  A  the  directions-of-arrival  of  the  signals  arriving 
at  the  array  need  to  be  estimated.  At  the  expense  of  more  com¬ 
putational  processing,  a  high-resolution  estimation  of  the  domi¬ 
nant  directions  can  be  obtained  applying  MUSIC  [21] ;  with  con¬ 
siderably  less  computation,  a  lower  resolution  estimation  can  be 
obtained  with  the  CBF. 

Substituting  A,  P^,  into  (35)  yields  a  parametric  estimate 
of  the  covariance  matrix  denoted  Ridmr-  To  avoid  the  look-di¬ 
rection  mismatch  a  signal  in  the  vicinity  of  the  look  direction 
should  not  be  included  in  the  formation  of  Ridmr-  The  MVDR 
beamformer  for  the  signal  si(t)  is  then  formed  using  Ridmr  as 


W idmr  — 


idmr 

al"Rirfmral  ' 


(38) 


TABLE  II 

Directions  and  SNRs  of  Incident  Signals 


Angle  (dg) 

-52 

-35 

0 

7 

12 

48 

SNR  (dB) 

16.2 

16.2 

-3.8 

16.2 

16.2 

16.2 

The  underlying  idea  of  the  indirect-DMR  technique  proposed 
herein  is  to  form  a  parametric  estimate  of  the  covariance  matrix 
given  information  on  the  location  and  power  of  the  dominant 
interferers  previously  obtained  from  either  MUSIC  or  the  CBF. 
With  this  parametric  estimate  of  the  covariance  matrix  where  the 
cross  correlations  are  removed,  the  MVDR  beamformer  maxi¬ 
mizes  the  output  SINR,  by  effectively  allocating  nulls  at  the  po¬ 
sitions  of  the  dominant  interferers.  It  is  important  to  note  that 
low  power  interferers  are  not  important  when  forming  the  para¬ 
metric  covariance  matrix  Ridmr,  since  it  is  not  important  to  form 
nulls  toward  them.  Thus,  the  performance  of  the  IDMR  beam- 
former  is  not  compromised  if  weak  sources  are  not  detected  by 
the  MUSIC  algorithm. 

There  are  several  works  that  focus  on  adaptive  beamforming 
with  low  sample  support.  Diagonal  loading  of  the  sample  covari¬ 
ance  matrix  has  been  extensively  used  in  adaptive  beamforming, 
however  there  has  always  been  the  issue  of  determining  the  op¬ 
timal  load.  Mestre  and  Lagunas  in  [13]  proposed  a  technique 
that  optimizes  diagonal  loading.  In  the  IDMR  algorithm,  there 
is  also  some  reduced  diagonal  loading  in  the  parametrically 
formed  covariance  matrix  Ridmr  (see  Table  I).  Vorobyov  et  al. 
proposed  in  [12]  a  robust  adaptive  beamforming  algorithm  using 
worst  case  performance  optimization.  In  Section  VI,  we  com¬ 
pared  the  performance  of  IDMR  to  the  performance  obtained 
with  Vorobyov’s  beamformer. 

VI.  Simulation  Results 

Simulations  were  conducted  employing  a  uniform  and 
linear  array  of  m  =  24  elements  with  half-wavelength  spacing 
receiving  plane- wave  signals.  In  the  simulations,  perfect  knowl¬ 
edge  of  the  array  manifold  was  assumed,  that  is,  no  steering 
vector  mismatch  was  considered.  It  is  relevant  to  note  that 
perfect  knowledge  of  the  array  manifold  does  not  imply  perfect 
knowledge  of  the  matrix  A.  To  estimate  A,  it  is  also  necessary 
to  know  the  directions-of-arrival  (DOAs)  of  the  interference 
signals,  and  in  the  simulations,  the  interference  DOAs  are 
estimated.  The  noise  at  each  array  element  is  spatially  and 
temporally  white  Gaussian.  The  incident  signals  are  modeled 
as  narrowband  with  amplitudes  modeled  as  complex  Gaussian 
random  processes.  We  considered  a  scenario  with  six  uncorre¬ 
lated  incident  signals  with  arrival  directions  (in  degrees)  and 
respective  signal-to-noise  ratios  (SNRs)  (in  decibels)  at  each 
array  element  as  shown  in  Table  II. 

Fig.  1  shows  the  output  SINR  versus  the  rank  of  the  beam- 
former  for  the  signals  at  0°  and  7°;  the  SINR  plotted  is  the  av¬ 
erage  over  200  Monte  Carlo  simulation  runs.  The  beamformers 
employed  were  the  SI-CG,  SD-CG,  and  DMR.  The  DMR  beam- 
former  requires  an  estimate  of  the  array  noise  power,  and  since 
this  beamformer  is  included  with  the  only  purpose  to  serve  as 
a  reference,  the  true  value  of  the  noise  power  is  used  as  its 
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(a) 


(a) 


Fig.  1.  SINR  performance  for  the  CG  and  DMR  beamformer  using  R  and  Fi§-  2-  SINR  Performance  for  low  power  signal  at  0°.  (a)  Six  snapshots, 

and  24  snapshots,  (a)  Low  power  signal  at  0°.  (b)  High  power  signal  ^  Twenty-four  snapshots. 

at  7°. 


own  estimate.  In  accordance  with  Section  IV,  where  we  proved 
the  equality  of  the  SI-CG  and  SD-CG  beamformers,  the  output 
SINR  obtained  with  these  two  beamformers  is  shown  to  be  the 
same.  The  MVDR  beamformer  with  the  true  covariance  matrix 
is  plotted  to  serve  as  a  reference  of  optimality. 

In  Fig.  1,  we  also  show  the  case  that  the  signal-free  sample 
covariance  matrix  R*+n  is  available,  and  we  plotted  the  output 
SINR  employing  the  DMR  and  CG  beamformers  (from  now  on, 
we  only  refer  to  the  CG  beamformer,  since  SI-CG  and  SD-CG 
yield  the  same  beamformer).  The  desired  signal  is  not  present  in 
R^+n.  It  can  be  observed  that  the  output  SINR  is  substantially 
higher  when  using  R^+n  than  when  using  R.  The  degradation 
of  the  output  SINR  when  using  R  is  not  caused  by  the  pres¬ 
ence  of  the  desired  signal  by  itself,  but  by  the  effect  of  residual 
correlations  between  the  desired  signal  and  the  interferers.  We 
observe  from  (34)  that  the  variance  of  the  residual  correlation 
is  proportional  to  the  power  of  the  desired  signal,  to  the  power 


of  the  interfere^  and  to  the  inverse  of  the  number  of  snapshots 
(1/n).  In  Fig.  1(a),  the  desired  signal  (located  at  0°)  has  a  power 
of  20  dB  (100  times)  less  than  the  desired  signal  (located  at  7°) 
in  Fig.  1(b).  We  observe  that  for  the  same  number  of  24  snap¬ 
shots,  the  degradation  when  using  R  instead  of  R^+n  is  much 
severe  in  the  case  of  Fig.  1(b)  (approximately  25  dB)  than  in  the 
case  of  Fig.  1(a)  (approximately  8  dB).  This  is  because  the  vari¬ 
ance  of  the  residual  correlations  in  Fig.  1(b)  is  100  times  larger 
than  in  Fig.  1(a). 

Figs.  2  and  3  show  the  output  SINR  when  the  IDMR  beam- 
former  is  used.  The  IDMR  beamformer  is  implemented  with 
both  the  CBF  and  MUSIC  to  locate  the  dominant  peaks.  While 
MUSIC  yields  finer  estimates  than  the  CBF,  the  former  requires 
significantly  more  computational  power.  We  observe  that  even 
when  the  CBF  is  used,  the  results  are  significantly  better  than 
those  obtained  using  the  existing  low-rank  beamforming  tech¬ 
niques,  such  as  CG  and  DMR.  The  reason  for  this  improvement 
is  that,  in  the  IDMR  technique,  the  residual  correlations  are  not 
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(a) 


(b) 


Fig.  3.  SINR  performance  for  low  power  signal  at  7°.  (a)  Six  snapshots, 
(b)  Twenty-four  snapshots. 


present  in  the  parametrically  formed  covariance  matrix.  It  can 
be  observed  that  the  residual  correlations  is  a  big  problem,  espe¬ 
cially  under  low  sample  support  (see  (34)).  We  observe  in  Figs.  2 
and  3  that  when  the  signal-free  covariance  matrix  R^+n  is  used, 
the  optimum  SINR  achieved  using  CG  is  close  to  the  SINR  ob¬ 
tained  with  the  IDMR  technique.  Under  higher  sample  supports, 
the  effect  of  residual  correlations  is  not  so  dramatic,  and  the  im¬ 
provement  obtained  with  IDMR  is  not  so  expressive.  The  rank 
of  the  IDMR  beamformer  is  the  number  of  estimated  signal  di¬ 
rections  used  to  form  the  matrix  A.  The  directions  are  sorted  by 
descending  order  of  estimated  signal  power;  that  is,  the  steering 
vector  of  the  signal  with  the  strongest  power  is  used  as  the  first 
column  of  A  and  so  forth.  The  rank  of  the  beamformer  is  equiv¬ 
alent  to  the  number  of  columns  of  A.  The  number  of  spectrum 
peaks  obtained  with  either  MUSIC  or  the  CBF  is  always  smaller 
than  M  =  24,  so  the  maximum  rank  of  the  IDMR  beamformer 
is  less  than  24.  The  figures  show  the  output  SINR  until  the  rank 


of  24.  Thus,  we  repeat  for  the  remaining  ranks,  the  output  SINR 
achieved  when  all  the  distinct  peaks  of  the  MUSIC  or  CBF  spec¬ 
trum  have  been  incorporated  into  A.  The  IDMR  beamformer 
requires  an  estimate  of  the  noise  power  at  the  array  elements;  in 
the  simulations,  we  used 

°n  =  A13  +  const  (39) 

where  A 13  is  the  thirteenth-largest  eigenvalue.  We  used  / 1000 
as  the  constant;  the  purpose  of  this  constant  is  to  not  let  <7^ 
vanish  when  the  number  of  snapshots  is  very  small.  When  using 
MUSIC  for  a  first  estimate  of  the  directions,  we  used  only  the 
12  noise  eigenvectors.  In  this  particular  case,  a  maximum  of  18 
could  have  been  used,  but  we  used  only  1 2  for  more  conservative 
results.  In  addition  to  a  higher  output  SINR  than  obtained  with 
either  CG  and  DMR,  IDMR  has  the  important  advantage  that 
the  maximum  SINR  over  all  the  ranks  is  very  close  to  the  SINR 
obtained  at  the  maximum  rank;  thus,  it  is  different  from  CG  and 
DMR  in  that  rank  selection  is  not  an  issue  (see  Figs.  2  and  3). 

Fig.  4  displays  the  plot  of  output  SINR  as  a  function  of  the 
number  of  snapshots.  For  the  CG  and  DMR  beamformers,  the 
plotted  SINR  value  is  the  average  over  200  Monte  Carlo  sim¬ 
ulation  runs  of  the  maximum  SINR  achieved  over  all  possible 
ranks.  For  the  IDMR  beamformer,  the  plotted  SINR  value  is  the 
average  of  the  SINR  at  the  maximum  rank.  It  is  important  to  note 
that  in  a  real-life  situation,  there  is  not  a  reliable  way  to  select  the 
optimal  rank  for  the  CG  and  DMR  beamformers.  Therefore,  the 
curves  show  an  upper-bound  performance  for  the  CG  and  DMR 
beamformers.  We  observe  that  the  IDMR  beamformer  yields  a 
high  SINR  even  when  a  very  small  number  of  snapshots  is  used, 
e.g.,  with  only  five  snapshots.  When  the  number  of  snapshots  is 
less  than  the  number  of  array  elements,  the  sample  covariance 
matrix  R  is  singular,  and  a  singular  value  decomposition  (SVD) 
was  used  to  estimate  the  noise  subspace. 

To  illustrate  the  effect  of  residual  correlations,  Fig.  5  shows 
a  plot  of  the  output  SINR  as  a  function  of  rank  for  the  case 
that  the  desired  signal  is  truly  correlated  to  an  interferer.  In  this 
scenario,  the  desired  signal  (located  at  0°)  is  partially  corre¬ 
lated  to  the  signal  located  at  7°  with  a  correlation  coefficient  of 
p  =  O.Se^/3.  Twenty-four  snapshots  were  used,  in  addition  to 
the  true  correlation  there  is  also  the  effect  of  the  residual  correla¬ 
tions  due  to  the  low  sample  support.  We  observe  that  even  when 
the  true  covariance  matrix  is  used  to  compute  the  MVDR  beam- 
former,  performance  is  still  degraded.  This  is  caused  by  the  true 
correlation  among  signals.  The  IDMR  beamformer  cancels  the 
correlation  embedded  in  covariance  matrix  and  the  output  SINR 
is  significantly  improved  even  when  the  CBF  beamformer  is 
used  to  estimate  the  directions  of  the  dominant  signals.  Observe 
that  when  using  the  CG  beamformer  with  the  desired  signal  free 
sample  covariance  matrix  R^+n,  the  output  SINR  is  similar  to 
that  obtained  with  the  IDMR  beamformer. 

In  Fig.  6,  we  compared  the  performance  of  the  IDMR  beam- 
former  and  the  robust  adaptive  beamforming  algorithm  (RAB A) 
proposed  by  Vorobyov  et  al.  in  [12].  RAB  A  is  known  to  be  ro¬ 
bust  against  low  sample  support  and  steering  vector  mismatch. 
To  facilitate  comparison,  we  used  the  same  scenario  which  was 
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.  4.  Output  SINR  versus  training  sample  size,  (a)  Signal  at  0° .  (b)  Signal  at 


Rank 

Fig.  5.  Output  SINR  for  truly  correlated  signals. 


Fig.  6.  Comparison  with  robust  adaptive  beamforming  algorithm  proposed  by 
Varobyov  et  al.  (a)  Without  look  direction  mismatch,  (b)  With  look  direction 
mismatch. 


used  in  [12].  We  considered  the  scenario  with  three  uncorre¬ 
lated  signals  arriving  at  the  array  with  directions  of  arrival  3°, 
30°,  and  50°  and  SNRs  of  — 10, 30,  and  30  dB,  respectively.  The 
array  is  linear  with  ten  omnidirecional  sensors  and  half- wave¬ 
length  spacing.  The  desired  signal  is  the  signal  arriving  at  3°. 
RABA  was  applied  with  the  constant  e  =  3.  In  Fig.  6(a),  there 
is  no  look-direction  mismatch,  and  the  output  SINR  is  plotted 
versus  the  number  of  snapshots.  In  Fig.  6(b),  there  is  a  look-di¬ 
rection  mismatch:  the  presumed  and  actual  directions  of  arrival 
are  3°  and  5°,  respectively.  In  this  simulation,  30  snapshots  were 
used  to  estimate  the  sample  covariance  matrix,  and  the  output 
SINR  is  plotted  versus  the  SNR  of  the  desired  signal  at  5° .  When 
MUSIC  is  used  in  the  first  step  of  IDMR,  the  output  SINR  ob¬ 
tained  with  IDMR  is  shown  to  be  higher  than  the  obtained  with 
RABA.  In  the  case  that  the  CBF  is  used  in  the  first  step  of  IDMR, 
then  RABA  yields  a  higher  SINR. 
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VII.  Conclusion 


Acknowledgment 


We  analyzed  existing  low-rank  beamformers  and  proved  that 
SI-CG  and  SD-CG  yield  the  same  low-rank  beamformer  at  any 
rank.  This  dictates  the  use  of  SI-CG  in  practice  since  it  does  not 
require  the  costly  computation  of  blocking  matrices.  Under  a  sit¬ 
uation  of  low  sample  support,  we  investigated  the  performance 
degradation  of  low-rank  MVDR  beamformers  due  to  residual 
correlations  embedded  in  the  sample  covariance  matrix.  This 
motivated  us  to  propose  a  new  adaptive  beamforming  technique 
(IDMR)  based  on  a  parametric  estimate  of  the  covariance  matrix 
in  which  correlations  among  signals  are  effectively  removed. 
Simulations  revealed  that  the  IDMR  beamformer  yields  a  dra¬ 
matic  improvement  in  the  output  SINR  relative  to  existing  CG 
and  PCI/DMR  beamformers,  and  unlike  the  CG  and  PCI/DMR 
beamformers  the  performance  obtained  with  the  IDMR  beam- 
former  does  not  depend  on  rank  selection. 

Appendix 

Blocking  Matrix  Effected  Through 
Householder  Transformation 

A  Householder  transformation  is  of  the  form 

H  =  Im  —  2vvff,  where  vHv  =  1.  (40) 

Im  is  the  to  x  to  identity  matrix.  It  can  be  shown  that  the  matrix 
H  is  an  unitary  matrix  for  any  unit  norm  vector  v.  Without  loss 
of  generality,  assuming  that  the  steering  vector  ai  has  a  unit 
norm,  we  consider  a  unitary  matrix  constructed  as 

H  =  [ai  !  B]  (41) 

where  B  is  referred  to  as  the  blocking  matrix,  with  B^ai  =  0 
and  B^B  =  Im_i.  It  follows  from  these  stipulations  that  H 
is  a  unitary  matrix  whose  “action”  on  the  steering  vector  is  as 
follows: 


H^ai  =6 1  = 


Substituting  (40)  into  (42)  implies  that 
vnai  — 


(42) 


(43) 


where  a  indicates  proportionality.  Defining  an  (m  —  1)  x  m 
selection  matrix,  F,  as 


r  =  [o ;  im_i]  (44) 

we  observe  that 

Bh  =  mH.  (45) 

Thus,  through  a  Householder  transformation,  the  blocking  ma¬ 
trix  B  can  be  computed  by 

c 

B  =  (Im  —  2vv^)Ft,  where  v  =  tt— ^ T.  (46) 

||ai  —  di|| 
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